clear;
load channel_task_4.mat

N=64;
M=16;

TrainingSize = 2;
MessageSize = 300;

FFT_NUM = 1024;
fn = 20e6*linspace(0,1,FFT_NUM);

hk1 = fft(h11, 64);
hk2 = fft(h12, 64);

mapping = @(input) (input * -2 + 1);
demapping = @(input) (real(input) < 0);

errors = 0;
total = 0;
for it=0:0
    messagesize = N*MessageSize;
    trainingsizebits = N*TrainingSize;
    data = randi([0, 1], 1, messagesize);
    %training = randi([0, 1], 1, trainingsizebits);
    %training = ones(1, trainingsizebits);
    r = [1 1 1 -1];
    r = [r r r -r];
    r = [r r r -r];
    training = repmat(r, 1, TrainingSize);
    
    b = [demapping(training) data];
    
    xk = mapping(b);
    [xm,xm2] = modulation(xk, M, N);
    
    ym1 = filter(h11, [1], xm);
    ym2 = filter(h12, [1], xm);
    
    %uncomment to remove channel filtering
    %ym1 = xm;
    %ym2 = xm;
    
    interference = pbj(size(ym1, 2), -10);
    interference1 = filter(hj1, [1], interference);
    interference2 = filter(hj2, [1], interference);
    
    %comment out to remove interference from bluetooth
    ym1 = ym1 + interference1;
    ym2 = ym2 + interference2;
    
    %comment out to remove awgn noise
    ym1 = awgn(ym1, 20);
    ym2 = awgn(ym2, 20);
    
    [yk1,yk3] = demodulation(ym1, M, N);
    [yk2,yk4] = demodulation(ym2, M, N);
    
    yout = zeros(1, messagesize);
    % apply beamformer algorithm to each subcarrier
    for k=0:N-1
        ryk1 = Ryk(yk1, yk2, k, TrainingSize, N)
        ryx1 = ryx(yk1, yk2, xk, k, TrainingSize, N);
        w = ryk1 \ ryx1;
        yout(k+1:N:end) = w' * [yk1(N*TrainingSize+1+k:N:end) ; yk2(N*TrainingSize+1+k:N:end)];
    end
    
    %uncomment to remove beamforming
    %yout = yk1(trainingsizebits+1:end);
    
    bout = demapping(yout);
    
    for k=0:N-1
        err = sum(abs((bout(k+1:N:end) - data(k+1:N:end))));
        fprintf('%d %d\n', k, err);
        %errors = errors + err;
    end
        
    errors = errors + sum(abs(bout - data));
    total = total + messagesize;
end

for i=2:10:100
    ryk1 = Ryk(yk1, yk2, 0, i, N);
    ryx1 = ryx(yk1, yk2, xk, 0, i, N);
    w = ryk1 \ ryx1;
end

BER = errors/total

%Pxx = abs(fft(ym,FFT_NUM)).^2/length(ym)/20e6;
%Hpsd = dspdata.psd(Pxx(1:length(Pxx)/2),'Fs',20e6);
%plot(Hpsd);
